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ABSTRACT 

A  nonlinear  thin  shell  theory  is  derived  for  the 
axl symmetric  buckling  of  spherical  shells  subjected  to 
either  a  pressure  or  a  centrally  directed  surface  load. 
The  theory  is  reduced  to  a  boundary  value  problem  for  a 
system  of  four  first  order  ordinary  differential  equations, 
Numerical  solutions  of  this  boundary  value  problem  are 
obtained  by  the  shooting  and  parallel  shooting  methods. 
An  extensive  numerical  study  is  made  of  the  nonlinear 
deformations  of  the  shells.   We  find  for  example,  that  all 
solution  branches  that  bifurcate  from  the  eigenvalues  of 
the  linearized  buckling  theory  are  connected  to  each  other 
by  means  of  intermediate  branches.   Some  Implications  of 
the  numerical  results  concerning  the  buckling  of  spherical 
shells  are  discussed. 
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1.   Introduction. 

The  outer  surface  of  a  thin-walled  hollow  sphere*  is 
subjected  to  a  uniform  compressive  load  p.   It  may  be  either  a 
pressure  or  a  centrally  directed  load.   The  tiniformly  contracted 
spherical  shell  is  a  possible  equilibrium  state  for  all  values 
of  p.   This  state  occurs  in  experiments  only  for  sufficiently 
small  p.   When  p  increases  to  or  beyond  some  critical  pressure 
the  shell  suddenly  Jumps  or  buckles  with  finite  displacements 
into  a  non-spherical  shape. 

The  classical  or  linear  buckling  theory  [1]  yields  an 
eigenvalue  problem  with  the  eigenvalue  parameter  proportional 
to  p.   The  spectrum  is  discrete  and  the  lowest  eigenvalue  is 
usually  called  the  buckling  load.   For  sufficiently  thin  shells 
the  eigenfunction  or  buckling  mode  consists  of  many  waves  over 
the  surface  of  the  shell  with  an  undetermined  amplitude.   The 
experimentally  measured  buckling  loads  are  frequently  consider- 
ably less  than  the  lowest  eigenvalue.   The  experimentally 
observed  deformation  is  usually  a  "dimple"  confined  to  a 
neighborhood  of  one  point  with  the  remainder  of  the  shell  remain- 
ing nearly  spherical.   Karman  and  Tsien  [2]  suggested  that  a 
nonlinear  theory  is  required  to  adequately  describe  the  buckling 
of  spherical  shells.   They  treated  the  complete  spherical  shell 
as  a  clamped  shallow  spherical  cap  corresponding  to  the  region 
where  the  dimple  occurs.   Their  analysis  was  subsequently  refined 
and  generalized  in  [3,^,5,6,7]  and  elsewhere.    The  solutions 


1 

For  example,  a  ping-pong  ball. 

We  refer  to  the  discussion  and  the  extensive  bibliographies 
in  [7j8]  for  references  to  previous  investigations. 
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of  a  nonlinear  theory  have  been  analyzed  in  [8,9,15]  t>y  a  perturba- 
tion methods  for  loads  near  the  eigenvalues  of  the  linearized 
theory . 

In  this  paper  we  consider  the  nonlinear  axisymmetric 
deformations  of  the  hollow  sphere.   An  appropriate  shell  theory 
is  derived  in  Section  2  and  in  Appendix  I.   The  strain  displace- 
ment relations  and  the  curvature  expressions  differ  from  those 
employed  in  [8,9].   The  shell  theory  is  reduced  to  a  two  point 
boundary  value  problem,  which  we  call  Problem  S,  for  a  system 
of  four  first  order  ordinary  differential  equations.   The 
boundary  value  problem  is  the  same  for  pressure  or  centrally 
directed  loading.   We  employ  the  shooting  and  parallel  shooting 
methods  to  obtain  extensive  numerical  solutions  of  Problem  S. 
The  numerical  methods  are  described  in  Section  4.   The  results 
of  the  computations  are  presented  in  Section  5.   Their  implica- 
tions concerning  buckling  are  discussed  in  Section  6. 

2 .   Formulation. 

Let  r ,  9 ,  (})  denote  a  spherical  coordinate  system. 
The  radial  distance  from  the  origin  is  r.   The  polar  angle  9 
is  measured  from  the  north  pole,  9=0,  and  is  in  the  range 
0  <_  9  <_  TT.   The  circumferential  angle  cf)  varies  in  the  interval 
0  <_<()<_  2 TT .   The  inner  and  outer  surfaces  of  the  thin  spherical 
shell  are  given  by  r  =  R  +  h  where  R  is  the  radius  of  the 
midsurface  of  the  shell  and  2h  is  the  uniform  thickness.   The 
shell  is  deformed  by  a  uniform  compressive  surface  load  p. 
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It  is  either  directed  toward  the  center  r  =  0  or  it  is 

a  pressure  that  is  normal  to  the  deformed  surface.   We  show 

in  Appendix  I  that  our  shell  theory  is  the  same  for  both  loads. 

We  consider  only  the  axisymmetric  deformations  of  the 
shell.   Thus  all  displacements  and  stresses  are  independent 
of  (() .   The  displacements  in  the  8  and  the  negative  r  directions 
are  denoted  respectively  by  U(r,0)  and  W(r,9).   The  (p   displace- 
ment vanishes  because  the  deformation  is  axisymmetric.   The 
exact  Lagranglan  strain  displacement  relations  in  spherical 
coordinates  are 


(2.1a) 


*         (U  -W) 
e««Cr,e)  =    \ 


1.1 


(U^-W)i 


+ 


■u+w,  ^ 


(2.1b)       eL(r.e)  =  '"  °°^^-  "' 


1   + 


1  (U  cot  9-  W) 


■]. 


(2.1c) 


e*e(r,0)  =  I 


r(U+W.)  _      _     U.-W 
^  (W  -1)  +  U  (1  +  -^- — ) 


Here  the  subscripts  on  e   signify  the  components  of  strain. 

The  subscripts  on  U  and  W  denote  partial  differentiation. 

The  shear  strains  e  ,  and  e„,  vanish  because  the  deformations 

rep       6<p 

* 

are  axisymmetric.   The  strain  e   is  not  included  since  it  is 

subsequently  neglected  in  deriving  the  shell  theory. 

We  assume  that  the  displacements  are  small  but  finite. 
Specifically,  this  means  that 


(2.2) 


Ug-w|    |u  cot  e-  w] 


w  I  <<  1 
r ' 
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We  use  (2.2)  in  (2.1)  and  obtain  the  approximate  nonlinear 
strain  displacement  relations. 


(2.3a) 


„  U„-W    ,  f\]+   W„^2 


(2.3b) 


*   0^  TT  ^   Q\  -  U  cot 
^AA  '^  E,(r,e)  =  : 


6-  W 


(2.3c) 


're 


^  E 


(r,0)  =  i 


(U+WJ 


U   - 


The  customary  assumptions  of  shell  theory  are  that: 
the  shell  is  thin,  i.e.  h/r  <<  1;  normals  to  the  midsurface 
remain  normal  to  the  deformed  midsurface;  the  normal  stress 
in  the  r  direction  is  negligible  compared  to  the  other 
normal  stresses;  the  absolute  values  of  the  strains  <<  1; 
the  two  dimensional  Hooke's  law  is  valid.  Thus  we  define  a 
new  dependent  variable  z  by 


(2.4) 


z  =  R 


z|  <  h  , 


and  we  assume  that 

(2.5)    W(r,e)  =  RW(e)  ,      U(r,6)  =  RU(9)  -  zV(e)  . 

Here  U  and  W  are  dimensionless  midsurface  displacements  and 
V  is  proportional  to  the  angle  of  rotation  of  a  tangent  to  a 
meridian.   Inserting  (2.5)  in  (2.3c)  and  setting  E  „  =  0  on 
z  =  0  we  obtain 
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(2.6)  V  =  U  +  W  , 

where  a  prime  denotes  differentiation  with  respect  to  9. 

By  substituting  (2.4)-(2.6)  in  (2.3a,b)  we  obtain, 
since  | z/R |  <<  1, 


(2.7) 


EQ(z,e)  =  e^(Q)    -    ZKg(e)  , 
E,(z,e)  =  e.(0)  -  ZK.Ce)  . 

<p       <p      <i> 


Here  the  midsurface  strains,  e^  and  e,,  and  the  curvatures, 
K„  and  K,,  are  defined  by 

(2.8a)       Cq  E  U'  -  W  +  I  V^  ,   e   =  u  cot  e  -  W, 
iP-.8h)  Kq  =  V'/R  ,    K.    =    (V/R)cot  6  . 

The  relations  (2.8a)  have  been  previously  employed  by 
Frledrichs  [4]  but  the  curvature  expressions  (2.8b)  differ 
from  those  in  [4].   Kolter  [8],  Thompson  [9]  and  others  use 
(2.8)  with  V  replaced  by  W.   The  strain  displacement  equations 
previously  used  in  analyzing  spherical  caps,  see  e.g.  [10], 
are  obtained  from  (2.8)  by  replacing  V  by  W  and  cot  9  by  1/9. 

The  shell  equilibrium  equations  and  boundary  conditions 
are  derived  in  Appendix  I  using  the  principle  of  minimum 
potential  energy.   They  are  reduced  to  a  boundary  value  problem 
for  the  following  system  of  four  first  order  ordinary  differen- 
tial equations: 

(2.9a)  y'(0)  =  f(y,9)  ,      0  <  9  <  u  . 
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Here  y(6)  is  a  four  dimensional  vector  with  components  y.(9), 
1  =  1,2,3,^  and  f(y,9)  is  the  four  dimensional  vector  with 
components  f . ,  1  =  1,2,3,^  defined  by 

f^  =  -(l-v)y^  cot  9  +  y^  +  (k  cot^9-P)y^  +  y^Y^^cot  9  , 

fp  =  yn  , 

(2.9b)    ^  _   ^     2  2 

f^  =  ypCcot  e-v)  -  y^  cot  9  -  y^  -  y^Ccot  9)/2  , 

2 
f^  =  (1-v  )y-j_/k  -  vy^^  cot  9  . 

The  parameters  P  and  k  are  defined  by 

(2.9c)        P  E  2"^(p/E)(R/2h)  ,    k  =  3""^(h/R)^  , 

where  E  is  Young's  modulus  and  v  is  Poisson's  ratio.  We  refer 
to  P  as  the  load.  The  components  of  y  are  defined  in  terms  of 
physical  quantities  by 

(2.10)  y^  =  m(9),  ^^    =    q(9),   y^  =  t(9),   y^^  =  v(9)  =  V 

where  m,  q,  t  are  proportional,  respectively,  to  the  radial 
bending  moment,  the  transverse  shear,  the  circumferential 
membrane  stress;  see  equations  (A. 2),  (A. 15c)  and  (A.l8)  in 
Appendix  I. 

Since  the  stresses  are  bounded  and  the  polar  axis  is 
a  line  of  symmetry  we  have  the  boundary  conditions 

(2.11)  u  =  w'=q  =  0   at   9  =  0,tt 

Here  u(9)  e  U(9)  and  w  is  defined  by 

(2.12)  w(9)  =  W(9)  -  Wq  ,   Wq  =  Cl-v)P  . 


-6- 


The  quantity  W^  is  the  radial  displacement  of  the  uniformly 
contracted  or  unbuckled  state.   From  (2.6),  (2.10)  and  (2.12) 

(2.13)  V  s  u  +  w'  , 

so  that  (2.11)  and  (2.10)  imply  that 

(2.14)  y^  =  y^  =   0   ,      0   =   0,   TT   . 

We  refer  to  the  boundary  value  problem  (2.9)  and  (2.l4)  as 
Problem  S. 

If  y-j  and  y^  are  eliminated  from  (2.9)  we  obtain  an 
equivalent  system  of  two  second  order  equations  in  y^  =  q 
and  y^^  =  v: 

(2.15a)  Lq  +  vq  =  -  V  -  (v^/2)  cot  0  , 

(2.15b)  Lv  -  vv  =  [ (l-v^)/k][-Pv  +  q  +  qv  cot  6]  , 

where  the  linear  differential  operator  L  is  defined  by 

(2.15c)  Lf  E  f"  +  f  cot  9  -  f  cot^e  . 

The  equations  of  the  classical  linear  theory  of  bending  of 
spherical  shells  [11]  are  obtained  from  (2.15)  by   setting  P  =  0 
and  omitting  the  nonlinear  terms  in  (2.15a,b).   It  is  not 
difficult  to  show  that  (2.15)  or  (2.9)  reduce  to  the  differential 
equations  of  the  nonlinear  spherical  membrane  [12]  when  we 
formally  set  k  =  0.   After  multiplying  (2.15b)  by  k  and  setting 
k  =  0,  the  left  hand  side  of  (2.15b)  formally  vanishes  if  9 
is  bounded  away  from  0  and  it  and  the  order  of  the  system  (2.15) 
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is  reduced.   The  reduction  is  not  uniform  in  6  because  of  the 
singularity  of  cot  9  at  9  =  0,Tr.   Thus  9  =  0,tt  are,  in  some 
sense,  turning  points  of  (2.15)  and  boundary  layers  may  form 
at  9  =  0  or  0  =  TT   (or  both)  as  k  ^  0.   Friedrichs  [i|], 
Altshuler  [5]  and  to  some  extent  Gabril'iants  and  Feodos'ev  [6] 
attempted  to  analyze  this  boundary  layer. 

It  is  easy  to  show  that  if   |u|  <<  |w'|  and  9  <<  1  (so 
that  V  is  negligible  compared  to  9   ),  (2.15)  can  be  reduced 
to  the  spherical  cap  equations  [10]. 

If  y  (9)  is  a  solutiori  of  the  hemisphere  problem, 

(2.16a)  y°   =  f(y°,0)  ,       0  £  9  <  tt/2  , 

(2.16b)  y°  =  y°  =  0    ,   at   9  =  0  ,  ^/2  , 


then 

(2.17)  y(6)  = 


y°(e)     ,   0  <  9  <  tt/2  , 
y°(TT-9)   ,   tt/2  <  9  <  TT  , 


is  a  solution  of  Problem  S.   Such  a  solution,  which  is  symmetric 
with  respect  to  the  equatorial  plane,  is  referred  to  as  a 
symmetric  solution.   The  hemisphere  problem  (2.l6)  is  also 
physically  significant.   The  boundary  conditions  (2.l6b)  are 
consequences  of  the  physical  conditions  u(Tr/2)  =  w'(Tr/2)  = 
=  q(TT/2)  =  0.   They  imply  that  the  edge  of  the  hemisphere 
cannot  rotate  or  deflect  normal  to  the  equatorial  plane  but 
it  can  slide  freely  tangential  to  the  plane. 

We  refer  to  solutions  of  Problem  S  which  do  not  satisfy 
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*       1 

y(S)  =  yCiT-O)  as  unsymmetric  solutions.    If  y  (6)  is  an 

un symmetric  solution  then  another  distinct  solution   is 

2       1  2  1 

y  (0)  =  y  (7T-0).   The  solution  y   Is  obtained  from  y   by 

rotating  the  sphere  so  that  the  north  and  south  poles  are 

Interchanged.   Thus  unsymmetric  solutions  occur  in  pairs. 

When  a  solution  of  Problem  S  has  been  obtained  u  is 

determined  from, 

(2.18)  u'  -  u  cot  9  =  (l+v)(q  cot  9  -  t)  -  v^/2  , 

and  w  from  (2.13).   Equation  (2.18)  is  singular  at  0  =  0,Tr. 
Therefore  the  boundary  conditions  u(0)  =  0    u('fT)  =  0  are 
insufficient  to  uniquely  solve  (2.l8).   Unique  solutions  of 
(2.18)  are  obtained  by  specifying  u'(0)  such  that  u(7r/2)  =  0. 
This  is  the  correct  boundary  condition  for  the  hemisphere  problem. 
It  fixes  the  rigid  body  displacement  in  the  complete  sphere  problem. 

The  difference  between  the  potential  energy  of  a  buckled 
solution  and  the  unbuckled  solution  is  proportional  to 


IT 

(2.19)  e  = 

0 


[t^  +  q^(cot^9-v)  +  (l-v^)m^/k  +  kv^cot^O  -4Pw]sin  0  d9 


* 

They  are  unsymmetric  with  respect  to  the  equator,  but  still 

retain  their  axial  symmetry. 
See  Appendix  I. 
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3 .   The  Linear  Buckling  Theory. 

The  linear  buckling  theory  is  obtained  by  omitting 
the  nonlinear  terms  (v  /2)  cot  6  in  (2.15a)  and  qv  cot  0  in 
(2.15b).   The  resulting  eigenvalue  problem  has  non-zero  solutions 
if  and  only  if 

(3.1a)   P  =  P^(k)  =  [k/(l-v^)][A^  +  1  +  v]  +  1/[A^  +  1  -  v], 

n  =  1,2,...  , 
where 

(3.1b)  X  =    n(n+l)  -  2  . 

The  corresponding  eigenf unctions  are 

(3.2)   q(0)  =  a^A^^^e)  ,   v(e)  =  B^A^^^G)  , 

^n  -  ^n^^^n'*'^"''^  '     ^  =  1>2,...  , 

where  3   t^  0  are  arbitrary  constants  and  A^   (G)  is  the  Associated 
Legendre  Function  of  the  first  kind  of  degree  n  and  order  1,  [13].  II 
n  is  even  the  corresponding  elgenfunction  is  symmetric  with  respect 
to  6  =  tt/2   and  if  P   is  simple,  i.e.  P   7^  P   for  any  n  7^  m,  it  is 
called  a  symmetric  eigenvalue.   If  n  is  odd  then  the  elgenfunction 
is  anti-symmetric  and  if  P   is  simple  it  is  called  an  unsymmetric 
eigenvalue.   The  value  of  n  for  which  P   is  a  minimum  depends  on 
k.   It  may  correspond  to  either  a  symmetric  or  an  unsymmetric 
eigenvalue. 

The  eigenvalues  obtained  from  previous  linearized  theories 
[1],  differ  slightly  from  (1).   For  example,  the  eigenvalues 
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determined  by  Zoelly  and  Schwerln  [1]  which  we  denote  by  P 
are  given  by 

2.  ^,  .  /.  .  .v2- 


(3.3)   P 


f 


IX    +2X    +(l+v)  ] 


n   n 


""  (X^+l+3v) 


k 


l-v2 


+  (X^+l+3v)"^,  n=l,2,. 


For  n  large  we  obtain  from  (3.1)  and  (3.3)  that 

^n=  ^nr^-^f^  (1^^)  7^+  ••• 
1-v     n         l~v 

2 

^n  =  ^n  ;^  -^  —  -^  [(1-3V)  +  ^-y-]  —^   +  ...  . 
1-v     n  n   1-v 

f 

Thus,  for  large  n,  P   is  only  slightly  larger  than  P   since 

k  is  small. 

In  the  numerical  computations  described  in  the  following 

sections  we  obtain  solutions  for  a  wide  range  of  P  values  for 

-3  -S 

k  =  10   and  k  =  10   .   In  Table  I  we  list,  in  increasing 

magnitude,  the  fifteen  lowest  eigenvalues  computed  from  (3.1) 

and  (3.3)  for  these  values  of  k. 

For  each  of  these  values  of  k  the  lowest  eigenvalue 

IS  unsymmetric.   There  are  other  values  of  k,  e.g.  k  =  1.2x10  "^ , 


for  which  the  lowest  eigenvalue  is  symmetric.   There  are  even 

^_  2 

values  of  k  =  k  _  E  ^^r — -r. — -y-T^r — -Tj — -y   for  which  the  corresponding 
mn    (A  +1-v)(a  +1-v)  '^  ^ 

n       m     ^ 

eigenvalue  is  a  multiple  eigenvalue.   For  example,  with 

k  =  .76935><10   ,  the  lowest  eigenvalue  occurs  for  n  =  5  and  n  =  6, 
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Table   lA ,    k  =   10 


-3 


Table   IB,    k  =   10 


-5 


n 

P 
n 

^n 

5 

.0675324 

.0646786 

6 

.0706160 

.0684843 

4 

.0750573 

.0702911 

7 

.0799192 

.0781080 

8 

.0936046 

.0919527 

3 

.106244 

.0949524 

9 

.110786 

.109220 

10 

.130993 

.129476 

11 

.153954 

.152467 

12 

.179504 

.178037 

13 

.207540 

.206085 

2 

.219602 

.172597 

14 

.237992 

.236545 

15 

.270812 

.269371 

16 

.305967 

.304531 

n 

P  xlOO 
n 

P  xioo 
n 

17 

.668365 

.665572 

16 

.671713 

.668555 

18 

.673789 

.671270 

15 

.685593 

.681940 

19 

.686669 

.684358 

20 

.706003 

.703853 

14 

.712402 

.708064 

21 

.731019 

.728995 

1 

|13 

.755470 

.750161 

22 

i 

.761113 

.759190 

23 

.795809 

.793966 

12 

.819535 

.812815 

24 

.834726 

.832946 

25 

.877555 

.875828 

11 

t 

.911529 

.902695 
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4 .   Numerical  Methods. 

For  both  the  sphere  and  the  hemisphere  problems,  we 
must  solve  (2.9)  subject  to  boundary  conditions  of  the  form 

(4.1)        y2(0)  =  yi|(0)  =  0  ,    y2(L)  =  y^(L)  =  0  . 

For  the  hemisphere  L  =  tt/2  and  for  the  full  sphere  L  =  it  . 
These  problems  are  solved  numerically  by  the  shooting  method 
for  the  hemisphere,  and  a  parallel  shooting  method  for  the  full 
sphere.   A  crucial  step  is  a  continuation  procedure  used  to  get 
good  estimates  of  the  proper  initial  data.   First  we  describe 
our  method  for  the  hemisphere  problem  and  then  indicate  how  that 
procedure  is  employed  to  solve  the  full  sphere  case.   The 
numerical  procedures  and  the  continuation  scheme  are  analyzed 
rather  thoroughly  in  [l4]. 

4 . 1   Shooting  for  Hemispheres. 

Consider  the  Initial  value  problem,  with  fixed  P  and  k, 

(4. 2a)  u'(e)  =  f (u,e)  , 


0 
(4.2b)  u(0)  = 


We  denote  the  solution  of  (4.2)  by  u  =  u(.Q  ;£,-.  ,E,    )  .      If  this 
solution  exists  in  0  <_  G  <_  tt/2  for  initial  parameters  ?,  and  5p 
such  that 
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(4.3)         \x^{ti/2;^^,^^)    =    0  ,   u^  (tt/2  ;^-^  ,^2^  =  0  ' 

then  by  (4.1),  y(9)  =  u{Q  ;E,      E,„)    is  a  solution  of  the  hemi- 
sphere problem.   The  number  of  distinct  solutions  of  the 
hemisphere  problem  is  equal  to  the  number  of  distinct  real 
roots  (C-|,Cp)  of  (4.3).   These  roots  are  numerically  determined 
using  Newton's  method,  and  a  continuity  procedure.   For  this 
purpose  we  introduce  the  vectors 

(4.4a)       U^^^e)  =  -^ri-   ,  u^^^e) 


Then  it  follows  by  differentiation  in  (4.2)  that  the  U^^ ^ 
satisfy  the  linear  variational  systems: 

dU^^) 
(4.4b)       -|g—  =  A(u,e)U^^\   U^^^O)  =  e^^^  = 


0 


dU^^^  (2)     (2)  (2)        '  °   ■' 

(4.4c)        -^.^  =  A(u,e)U^'^\   U^^^O)  =  e^'^''  =  !      i 


de 


0 


Here  A  is  the  fourth  order  Jacobian  matrix, 

/Sf. (u,e) 

C4.4d)  A(u,e)  =      -^^ 

~       \     o  u . 

Given  an  initial  estimate  (c|   ,5p  ^)    for  a  root  of 
(4.3),  the  sequence  of  Newton  Iterates  {5n   jCp   }  is  defined  by 
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J  J  J 

where  (A^-^  ,    A?^   )  i^  ^he  solution  of  the  linear  system 


(4.5b) 


The  coefficient  elements  U.  '^(7t/2)  in  (4.5t))  are  the  corres- 
ponding  elements  in  the  solution  of  (4.4)  in  which 
u(9;C!   jCo   )  is  used  in  the  evaluation  of  A(u,9). 

A  Runge-Kutta  method  is  employed  to  solve  the  three 
initial  value  problems  (4.2),  (4.4b)  and  (4.4c).   The 
singularity  at  9  =  0  requires  special  starting  procedures. 
In  (4.2a)  we  simply  use  the  limiting  relations 


(4.6) 


f(u(0),  0)  = 


0 

0 
1-v 


k    ''Ij 

obtained  from  (2.9b)  and  (4.2)  by  requiring  continuous 
derivatives  at  9  =  0 .  In  (4.4b)  and  (4.4c)  the  first. step  in 
the  integration  is  done  with  the  centered  scheme: 

(4.7)    U(h)  =  Tj(o)  +  ^  A(ku(0)  +  u(h)],  |)-[U(0)  +  U(h)] 

where  h  is  the  step  width.   The  error  is  0(h)  since  we  use 
this  second  order  scheme  only  at  one  point. 

The  interval  0  <  9  <  tt/2  is  divided  into  100  equal 
subintervals  and  the  integrations  are  performed,  starting  at 
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0  =  Qj  over  each  sublnterval  with  spacings  h„  =  7r/200,  and 
h-,  =  h„/2.   If  the  maximum  relative  difference  In  the  two 

Integrations  at  the  end  of  a  sublnterval  Is  not  less  than 

-M 
10   J  we  then  halve  the  step  width,  and  repeat  the  procedure 

over  the  sublnterval.   When  agreement  Is  achieved  we  proceed 

to  the  next  sublnterval,  starting  with  the  original  mesh 

widths  h„  and  h-,  .   The  Iterates  In  Newton's  method  (4.5)  are 

terminated  when  max  |A^.   /^.  '^1  <  10~  .   In  most  of  the 

1=1,2 
calculations  M  =  7-   Occasionally  a  larger  value  of  M  Is  used. 


4 . 2   Continuation  Method  for  Parameter  Variation. 

After  a  solution  Is  determined  for  some  value  of  P 
with  k  fixed,  (or  for  some  value  of  k  with  P  fixed)  accurate 
estimates  of  the  initial  Iterates  ?.  '^  to  use  for  P  +  AP  Cor 
for  k  +  Ak)  can  be  obtained  by  finding  9^./9P  (or  d^./dy:). 

J  J 

With  k  fixed  let  ?.(P)  denote  the  exact  initial  values  for 

J 

any  given  P.   Then  the  boundary  conditions  (4.3)  yield  the 
identities  in  P: 

u  (TT/2;C-,(P),ep(P),P)  -  0  , 
(4.3')        2      12 

Ui^(Tr/2;^^{?),i^{?),?)    =   0  . 

In  rewriting  (4.3)  we  have  Indicated  the  dependence  of  u  on  P 
which  was  previously  ignored  to  simplify  the  notation.   Differ- 
entiating each  of  these  identities  with  respect  to  P,  we  obtain, 
recalling  (4 .4a) , 
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V^2^h-n/2)k-^i?)    +   U^2^(TT/2)i2(P)  =  "  U^5)(^/2)  ^ 


(4.8) 


(1) 


(2) 


(5) 


UJ^^(tt/2)?^(P)  +  U^^^(tt/2)C2(P)  =  "  UJ^^(tt/2)  , 


where  we  have  Introduced 
(5) 


9u     .        3C.(P) 


(11. 9a)   U^^^e)  =  ^   ,   q(P)  =    ^p    ,   i  =  1,2. 

By  differentiating  in  (4.2)  we  find  with  the  aid  of  (2.9b) 
and  (4.4d)  that 


(4.9b) 


dU^^^e) 
de 


=  A(u,e)u*^5^  + 


^-^4(9)] 

0 

0 

0 

u'^5)(o)  =  0  . 


The  coefficients  in  (4.8)  are  already  computed  for 
use  in  the  final  Newton  iterate  (4.5b).   Thus  only  the 
additional  system  (4.9b)  is  solved  to  evaluate  ^^      and  ^^ 
from  (4.8).   Then  for  initial  estimates  we  use  the  first 
two  terms  in  the  Taylor  expansion: 


(4.10) 


k\^\f    +  AP)  =  C.(P)  +  ?.(P)AP 


J  =  1,2, 


The  derivatives  Ki'P)    are  of  physical  interest,  especially 

J 

when  P  is  near  a  bifurcation  point  or  near  a  local  maximum 
or  minimum  point  of  P(5.). 

The  corresponding  continuation  analysis  and  computa- 
tion for  variation  of  the  geometric  parameter  k  follow  in  an 
obvious  manner.   In  fact,  to  replace  P  by  k  only  the  right 
hand  side  in  (4.8)  need  by  changed  by  using  the  components  of 
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(4.11a) 


Proceeding  the  same  way  as  In  the  derivation  of  (4.9b)  we 
find  that 

U|^(e)cot  e 

0 


(4.11b) 


du<^^^e) 


de 


A(u,e)u^^^(e)+ 


0 


l-v' 


k 


2    1 


u,  (0) 


,  U^^^O) 


=  0. 


4.3   Parallel  Shooting  for  the  Full  Sphere. 

We  treat  the  full  sphere  by  employing  the  previous 
procedure  for  two  hemispheres.*   We  integrate  from  each  pole 
to  the  equator  and  adjust  four  parameters  in  order  to  satisfy 
four  continuity  conditions  at  the  equator. 

Consider  the  two  initial  value  problems. 


It  is  inconvenient  to  treat  the  full  sphere  by  Integrating  from 
9  =  0  and  adjusting  ^-,  and  Cp  to  satisfy  the  two  boundary 
conditions  at  6  =  tt.   Some  of  the  derivatives  are  indeterminate 
at  0  =  TT,  and  the  derivatives  of  some  of  the  auxiliary  variables 
U^    and  U^    are  infinite  there.   Numerically,  these  singulari- 
ties are  handled  more  readily  by  the  special  starting  procedure, 
(4.6)  and  (4.7),  than  by  integrating  up  to  ir-e  ,  and  then  using 
special  formulas. 
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(4.12a) 


(4.12b) 


du 

dt  =  ^^'i'^^    ' 


dv 

de 


u(0)  = 


f  ^1 

0 

?2 

0 

•                        • 

=  f(v,0)   ,     v(tt)  = 


0 
0 


,  0  <  e  £  tt/2; 


IT  >  0  >  it/2 


Here  the  E, .    are  the  four  initial  parameters,  two  for  each  of 
the  two  initial  value  problems.   The  solutions  of  these  problems 
are  denoted  by  u{e  ;^^,^^)    and  v{Q  •,^^,E,^)  .      If  we  can  find  C^ 
such  that 


(4.13)      Uj(tt/2;C^,C2)  =  Vj(tt/2;52,?^)  , 


J  =  1,2,3,4, 


then  a  solution  of  the  full  sphere  problem  is  given  by 


(4.14) 


y(e) 


_  j  u(e;C^,?2)  '     0  1  0  1  ^/2, 


For  fixed  P  and  k,  the  number  of  distinct  solutions  of 

Problem  S  is  equal  to  the  number  of  distinct  real  roots  {E,-.,E,2, 

?3,?^)  of  (4.13). 

The  equations  (4.13)  are  also  solved  by  Newton's  method. 
Now  however  there  are  two  nonlinear  systems  (4.12)  and  four 
linear  variational  systems  to  be  solved  for  each  Newton  itera- 
tion step.   The  numerical  integration  of  the  initial  value 
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problems  is  done  In  the  same  way  as  for  the  hemisphere.   All 
convergence  criteria  are  the  same  as  in  the  hemisphere  calcula- 
tions.  Corresponding  continuity  methods  in  both  k  and  P  are 
also  used  for  obtaining  initial  Iterates. 

4.4   Parallel  Shooting  for  Hemispheres. 

For  some  of  the  solution  branches  for  the  hemisphere  problem 
boundary  layers  occur  at  9  =  0  or  9  =  17/2;  see  the  discussion 
in  Section  5.   When  a  boundary  layer  occurs  at  9  =  tt/2  the 
solution  varies  rapidly  in  a  small  neighborhood  of  0  =  7t/2. 
Thus  small  variations  in  the  parameters  5-,  and  E,^   produce  large 
variations  in  the  solution  of  (4.2)  at  9  =  ■t\/2.      When  a  boundary 
layer  occurs  at  6  =  0,  Up  and  u^,  essentially  vanish  in  an 
interval  a  <  9  <_  tt/2.   Then  relatively  large  changes  in  E,^   and 
5p  produce  small  changes  in  the  solution  at  9  =  tt/2. 

In  both  cases  the  simple  shooting  method  of  Section  4.1 
is  unsatisfactory  for  an  accurate  determination  of  the  roots  of 
(4.3)  and  parallel  shooting  is  employed.   Thus  we  select  a  point 
Qr.   in  the  interval  0  <  9^^  <  tt/2  and  integrate  from  9  =  0  to  9^ 
and  from  9  =  tt/2  to  9^.   The  corresponding  four  initial  param- 
eters are  determined  from  continuity  requirements  analogous  to 
(4.13).   The  procedure  is  completely  analogous  to  the  parallel 
shooting  for  the  full  sphere  described  in  Section  4.3. 
Appropriate  values  of  9^,  are  determined  from  numerical  experi- 
ments . 
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5.   Presentation  of  Results. 

Extensive  numerical  solutions  of  Problem  S  have  been 

obtained  for  shells  with  k  =  10  ^  (the  thick  sphere)  and  k=10 

_5 
(the  thin  sphere)  over  a  wide  range  of  P  values.   For  k=1.2xl0 

a  less  extensive  P  range  has  been  covered.   At  the  fixed 

pressure  P  =  6xlO~   solutions  have  been  obtained  for  various  k 

in  the  range  10  -^  <_  k  <  2xl0~  .   We  shall  describe  the  thick 

sphere  results  that  seem  most  relevant  to  buckling  and  briefly 

mention  significant  differences  for  the  thin  sphere  case.   The 

implications  of  these  results  for  axlsymmetric  buckling  are 

discussed  in  Section  6. 

The  thick  sphere  computations  are  summarized  in 

Figures  1  and  2  where  we  plot  P  vs.  A(P)  and  e(P).   The  quantity 

IT 

2 

6 


(5.1)  A(P)  s  (i 


P     1/2 

w   d0)     >  0 


is  a  measure  of  the  amplitude  of  the  deviation  of  the  deforma- 
tion at  load  P  from  the  uniformly  compressed  state,  w  e  0. 
Similarly  e(P),  defined  in  (2.19),  is  proportional  to  the 
difference  between  the  potential  energies  of  the  buckled  and 
unbuckled  states  at  the  same  load,  P.   The  circles  and  crosses 
on  the  P-axis  denote,  respectively,  the  symmetric  and  unsymmetric 
eigenvalues  of  the  linear  theory  (see  Table  lA  in  Section  3). 

Each  point  (A(P),P)  with  A(P)  >  0  on  the  graphs  in 
Figure  1  represents  buckled  states,  i.e.  nontrivial  solutions 
of  Problem  S.   If  the  point  corresponds  to  a  symmetric  state. 
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then  it  represents  a  single  solution.   If  It  corresponds  to 

an  un symmetric  state,  then  It  represents  a  pair  of  solutions. 

The  graphs  show  that  for  fixed  values  of  P,  Problem  S  may  have 

many  solutions.   To  describe  and  discuss  this  multiplicity, 

we  decompose  the  graphs  of  M?)    Into  single  valued  segments 

of  maximal  extent  and  denote  each  such  segment  by  B   .  with 

n  J  J 

integers  n  and  j .   The  procedure  for  choosing  appropriate 

values  of  n  and  J  Is  described  in  Section  5.1.   This  is 

equivalent  to  decomposing  the  graph  of  e(P)  in  Figure  2  into 

monotone  segments  of  maximal  extent.   The  P  coordinates  of 

the  end  points  of  the  segments  are  denoted  by  P   .  and  P„  .  . 
^  n,j      n,j 

We  call  them,  respectively,  the  lower  and  upper  critical  points 
of  the  segment.   If  each  point  on  B   .  corresponds  to  a 

symmetric  [unsymmetrlc ]  solution,  then  we  denote  the  branch 

+ 

of  solutions  corresponding  to  E   .  as  B   .(P)  [B   .(P)]. 

^     ^     n,j      n,j  '      n,j 

The  P  axis  corresponds  to  the  unbuckled  states.   We  shall  not 
Include  this  branch  in  our  discussion  of  the  multiplicity. 

Buckled  solutions  of  Problem  S  are  said  to  bifurcate 
from  the  unbuckled  solution  at  or  from  a  pressure  P„  if  there 
is  a  solution  branch,  y(6,P)  depending  continuously  on  P  such 
that  yOjPg)  =    0  and  y(e,P)  ^  0  for  P  /^  Pg  t)ut  |P-Pq|  suffi- 
ciently small.   This  solution  branch  is  called  a  bifurcation 
branch.   The  bifurcation  branch  need  not  exist  for  P  in  a  full 
neighborhood  of  Pq.   If  It  exists,  locally,  for  P  >  Pq  (<  Pq) 
we  say  that  solutions  bifurcate  or  branch  up  (down)  from  P„. 

Our  calculations  show  that  there  are  solutions  In  a 
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full  neighborhood  of  each  of  the  symmetric  eigenvalues  that 
were  analyzed  and  that  solutions  branch  either  up  or  down  from 
each  of  the  unsymmetric  eigenvalues  that  were  analyzed.   Further- 
more the  results  indicate  that  ^   =  (dA(P  )/dP)~   =  0  for  n  odd, 
i.e.  unsymmetric  eigenvalues,  but  not  for  n  even  and  <^    is 
independent  of  k.   These  results  are  in  agreement  with  the 
conclusions  of  bifurcation  theory  [15]. 

5.1   The  Thick  Sphere,  k  =10"^. 

Perhaps  the  most  striking  observation  in  the  thick  sphere 

case  is  the  fact  that  all  computed  branches  of  solutions,  21 

symmetric  branches  and  16  pairs  of  unsymmetric  branches,  are 

connected.   Furthermore,  for  each  branch  B   .  the  corresponding 

n ,  J 

potential  energy  e   .(P)  is  a  monotone  decreasing  function  of  P 
on  P   •  1  P  1  P   •  (see  Figure  2a).   All  but  two  of  the  symmetric 

branches  and  five  pairs  of  the  unsymmetric  branches  contained 

F  E 

pressures  P  =  P   .,  called  equal  energy  loads,  at  which  e_  • (P^  .)=0 

for  a  buckled  state.   The  critical  pressures,  equal  energy  loads 

and  connectivity  of  the  unsymmetric  branches  and  the  first  16 

symmetric  branches  of  solutions  are  summarized  in  Table  II. 

Each  branch  bifurcating  from  an  eigenvalue  is  indicated  by  an 

asterisk.   This  table  complements  Figures  1  and  2. 

The  branch  B^  -,  bifurcating  up  and  down  from  the  lowest 

symmetric  eigenvalue  P^  =  0.7062  extends  from  the  lower  critical 

load  Pg  -L  -  .01045  to  the  upper  critical  load  Pg  -]_  -  .07127. 

At  its  lower  endpoint  B^  ^  Joins  the  branch  we  call  B^  „  whose 
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Table  IIA 


Symmetric  branches 


Branch  B 

pL        i 
n,j      ' 

P^  .  or  (P  ) 
n,,i     '  n' 

pU  . 
n,,l 

^6,0 

1 

,0104496    i 

.01214 

9 

^^6,1 

i 

(.0706150) 

.0712660 

^6,2 

.067830')    1 

.07123 

*   ^4,1 

! 
! 

C.O750573) 

.0855778 

^4,2 

! 
! 

.04^220^ 

.0848 

^4,3 

i 

"' 

.0662950 

^4,iJ 

i 



.04"^   22     i 

*   ^8,1 

(.0936046) 

.109261 

^8,2 

.1077 

.0656532 

^8,3 

.1358 

1 

.1366638 

*  ^10,1 

(.130993) 

.0735912 

*   ^12,1 

(.179504) 

.179807 

^12,2 

.105326 

.1798 

*  ^2,1 

(.219602) 

3 

.306228 

^2,2 

1  .237956 

.23797 

*  ^14,1 

1 
i 

(.2379"2) 

....  ij'^.r'  -.  ,, 
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Table   IIB 


Unsymmetric  Branches 

1            t        i 

'  Branch  B   .    ! 

^n,.1 

Pn,.1  °-  (Pn) 

n,.1 

1                               1 

i  *  B            'i 

,                                                                                    1 

.0128125 

(.0675324) 

.0675324 

"5.2 

.0137 

.0202242 

5,3        ! 

.0173^426 



■  ^,4        1 

{ 

.0199705 

\       "5,5        1 

.0199473 

1   "5,6 

.0208548 

.0123971 
(on  Bg^^) 

.01415 

■                   { 

.0799192 

(.0799192) 

1 

.0858161 

1   "7.2 

t 

.0467240 

;  .0855 

;'     ^7,3 

j 

i 

.064663 
(on  B^^^)  i 

'  "^3,1        i 

.106244 

— 
(.1062^4) 

1 

\ 

.110025   ; 

: 

1 

.105617 

.109531 

1 

IB,,         ; 
3,3        \ 

( 

.105790 
(on  Bg^2^ 

^^,1        ; 

.110786 

(.110786) 

> 
.120483   1 

^9,2 

i 

.0706469 

.118910 

^9,3 

1 

1 

1 

.10349    I 
(on  63^3)  ^ 
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lower  endpolnt  is  at  P^  „  =  P^  -,  and  whose  upper  endpolnt 

has  not  been  determined,  but  is  greater  than  P  =  .1;  see 

Figure  lb.   At  its  upper  endpolnt  B--  ,  joins  Bg  „  whose  upper 

endpolnt  is  at  P^  p  =  P.-  -,  and  whose  lower  endpolnt  is  at 

Pg  2  -  .06783;  see  Figures  la  and  2a.   The  labeling  of 

adjacent  joining  branches  follows  the  arbitrary  rule  implied 

above  with  the  proviso  that  whenver  a  bifurcation  branch  occurs, 

say  bifurcating  from  P  =  P  ,  we  call  it  B   ...   Thus  for  example, 

the  branch  joining  B^  „  at  its  lower  critical  pressure  should 

be  called  Br    o'   But  since  Br    n   bifurcates  from  P^  we  also  call 

it  Bk  -,  .   The  notation  in  Table  II  is  now  rather  obvious  if 

we  merely  observe  that  most  critical  pressures,  P   .  or  P   ., 

have  dual  representations.   On  the  branches  bifurcating  from 

unsymmetric  eigenvalues  the  eigenvalues  are  also  upper  or  lower 

critical  pressures,  that  is,  P   =  P   .  or  P   =  P   -,  . 
'^  '  '      n  n,l      n     n,l 

In  Figures  3  and  ^  we  show  how  the  solutions  for  P 
decreasing  on  the  lowest  symmetric  bifurcation  branch  Br    -,  go 
over  continuously  to  the  solutions  for  P  increasing  on  the 
adjacent  branch  B^  „.   There  is  an  obvious  difference  in  the 
shapes  and  stresses  between  the  states  on  branches  Br   -,  and  Br    q 
for  P  >  0.05;  compare  Figures  3a  and  ^a  with  Figures  3c  and  4c. 
However,  for  P  near  Pg  -,  -   ^ ({    n   ~    0.0105  the  shapes  on  both 
branches  are  quite  similar.   Thus  there  is  no  abrupt  change  in 
mode  (see  Figures  3b  and  4b).   As  P  Increases  the  states  on 
Br   Q   develop  a  boundary  layer  near  9  =  tt/2.   For  example,  at 
P  =  0.10,  t(0)  is  almost  constant  over  0  <  G  <  tt/3  and  then 
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changes  rapidly  over  it/3  <  9  <  77/2  to  satisfy  the  boundary 
conditions.   The  graph  of  w(0)  Is  approximately  a  parabola  with 
maximum  at  G  =  0  so  that  the  deformed  sphere  Is  essentially 
two  pushed  In   hemispheres  joined  at  the  equator.   The  develop- 
ment of  the  boundary  layer  and  the  simple  pushed-ln  shape  of 
the  sphere  on  B^  „  suggests  that  this  mode  and  branch  of 
solutions  continue   to  exist  as  P  ->■  +°°.   The  rapid  decay  of 
energy  e^  n^^^  would  seem  to  be  consistent  with  this  assumption. 

In  Figures  5  and  6  we  show  four  solutions  on  the 
branches  B^  .,  ,  B^  „  and  B^  ^  =  Bk  .,  .   One  of  these  solutions  Is 
for  P  =  .0712594,  near  the  critical  pressure  P^  ,  =  P^  pj  ^'^^ 
another  Is  for  P  =  .0678672  near  the   critical  pressure 

P^  n  =  P),  1  •   The  other  solutions  are  for  P  near  the  bifurcation 
6,2     4,1 

loads  ?r   and  ?u.      These  solutions  have  distinct  mode  shapes  and 

relatively  small  amplitudes  since  P  Is  close  to  P^  or  P/-.   This 

behavior  Is  typical  for  all  even  n  since  the  solutions  on 

different  symmetric  bifurcation  branches  B   .,  have  distinctly 

n ,  1 

different  deformation  and  stress  modes  for  P  near  P  ,  and  as  P 
varies  the  solutions  continuously  change  on  each  branch  and 
merge  with  the  solutions  on  some  other  branch  which  Is  usually 
not  a  bifurcation  branch.    Thus  all  symmetric  solutions  are 
connected  as  P  varies  up  and  down  on  the  appropriate  symmetric 
branches.   The  energy  variation  on  all  branches  Is,  as  previously 
noted,  monotone  decreasing.   However,  for  P  <  Pp  the  energy  on 


The  only  symmetric  bifurcation  branches  to  merge  directly  with 
each  other  were  B-,  „  .,  and  B,  „  ^  . 
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most  symmetric  branches  Is  mainly  positive  but  for  P  >  Pp  the 
opposite  occurs  (the  notable  exception  Is  B^  „,  see  Figure  2a 

and  Table  II) . 

+ 
The  pairs  of  unsymmetrlc  bifurcation  branches  B 

n,l 
merge  continuously  at  the  eigenvalues  P  =  P   for  n  odd.   The 

i  ^ 

branches  B^-  ,  bifurcate  down  from  the  lowest  eigenvalue 

±-       L 
Pf-  =  0.0675  and  merge,  respectively,  with  B^-  p  at  P^.  p  =  0.0128. 

Similarly,  Bj-  p  Joins  B^  ^   which  joins  B^-  k  which  joins  B^    ^ 

which  joins  B^    r   which  joins  B^-  y  which  finally  joins  Bj-  ^  at 

Pj-  y  =  0.0124;  see  Figure  2b  and  Table  IIB.  The  obvious  mergers 

also  occur  for  the  B^-  ..   Since  B^-  n   and  B~  „  loin  at  P;-  „  i.e., 

5, J         5,7      5,7  ^        5,7 

y'^(e,P^  ^)=  y"(Tr-e,P^  ^)  =  y^(TT-e,P^  ^),  B*  ^(P)  must  be 
symmetric  at  P  =  Pj-  „.   Furthermore,  the  results  show  that  this 
symmetric  state  lies  on  the  branch  B^  ^  .   The  value  P,-  „  is  not 


a 


critical  pressure  for  B^  .,  .   In  fact,  we  find  that 


P/-  -j  <  Pc-  y  <  Pr  -1  ;  see  Table  II.   Equal  energy  loads  are  found 

+        + 
on  the  branches  B^-  p  and  B^^  „.      Graphs  of  typical  solutions  on 

branch  B~  -,  are  shown  in  Figures  7a  and  8a.   Solutions  at  some 

of  the  upper  and  lower  critical  loads  P^-  .  and  P^  .  are  shown 

5, J      5, J 

In  Figures  7b  and  8b.   The  transition  from  one  branch  to  another 

Is  Indicated  by  these  graphs.   We  note  that  on  Bj-  ,  the  "southern" 

hemlsphere  is  less  distorted  than  the  "northern"  hemisphere  which 

suffers  a  relatively  large  indentation,  or  dimple,  around  the 

pole.   Similar  features  occur  in  all  the  B^  .  branches. 

5, J 

At  each  of  the  unsymmetrlc  eigenvalues  P  ,  for  n  =  7,3 

+ 

and  9  the  solutions  B   -,  bifurcate  up  and  go  through  a  series 

n ,  1 
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of  mergers  with  Intermediate  branches.   Finally,  B   .  Joins  B~  . 

for  some  j  at  P  =  P   .  or  P  =  P   . .   Again  the  symmetric  state 

Into  which  the  unsymmetrlc  pairs  of  solutions  finally  merge  Is 

found  to  lie  on  some  symmetric  branch  of  solutions.   For 

example,  B~  ^   join   ^u    h      at  P„  ^  =  0.0646.   Thus  we  find 

the  amazing  result  that  all  bifurcation  branches,  for  both 

symmetric  and  unsymmetrlc  modes,  are  connected  to  each  other 

by  means  of  Intermediate  branches.   The  energies  on  the 

bifurcation  branches  B~  ^  for  n  =  7,3,9  are  negative  for  P  >  P_. 

They  are  also  negative  on  B   „  for  a  limited  range  of  P  values. 

n ,  <i 


5.2   The  Thin  Spheres. 

_5 
The  results  for  the  thin  sphere,  k  =  10   ,  are  summarized 

In  Figures  9  and  10.   For  this  value  of  k  the  lowest  eigenvalue 

Is  also  unsymmetrlc.   The  eigenvalues  for  k  =  10     (see  Table  IB 

of  Section  3),  are  closely  spaced.   For  example  a  10^  change  In 

P  from  P  =  P-,y  spans  seven  additional  eigenvalues.   A  1%    change 

Includes  two  additional  eigenvalues.   The  solutions  for  the  thin 

sphere  are  qualitatively  similar,  in  many  respects,  to  those 

previously  discussed  for  the  thick  sphere.    All  the  branches 

bifurcating  from  symmetric  eigenvalues  are  continuously  connected 

in  the  same  way  as  for  the  thick  sphere.   The  solution  branch 


* 

Perhaps  the  most  significant  difference  is  the  Increase  in 

the  multiplicity  of  the  solutions  as  k  decreases. 
For  example,  for  the  thick  sphere  7  different  symmetric 
solutions  were  determined  at  P  =  .07.   For  the  thin  sphere 
46  different  symmetric  solutions  were  determined  at  P  =  .007 
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bifurcates  down  from  the  lowest  eigenvalue  Pny.   The  solution 
branches  bifurcate  up  from  the  next  two  unsymmetrlc  eigenvalues 
P-|  [-  and  P-,  Q  and  have  negative  energy  for  P   near  the  eigen- 
values.  We  have  not  completely  Investigated  the  unsymmetrlc 
branches  of  the  solution  for  the  thin  spheres.   However, 
numerical  results  that  are  not  shown  In  the  graphs  Indicate 
that  the  unsymmetrlc  solutions  connect  with  the  symmetric  ones 
In  the  same  way  as  for  the  thick  sphere.   The  formation  of 
boundary  layers  as  P  varies  on  B, ^  .  is  more  pronounced  for 
the  thin  sphere  than  for  the  corresponding  branch  of  the 
thick  sphere. 

A  feature  of  the  thin  sphere  that  has  not  been  observed 
for  the  thick  sphere  is  the  existence  of  Isolated  solution 
branches.   They  are  branches  that  are  not  connected  to  any  of 
the  solutions  which  bifurcate  from  the  eigenvalues.   For  example, 
one  of  these  occurs  for  P  in  the  interval  .OO6OO5  <  P  <  .010338. 
There  are  critical  loads   at  P  =  .007760,  .006095,  .010338 
and  .006005.   The  energies  on  this  branch  are  in  the  interval 
+  .119  <  (e/2)  ^ -^   <    .232.   Since  the  energy  is  large  and  positive 
this  branch  is  not  shown  in  Figure  9.   We  have  not  found  any 
Isolated  branches  with  negative  energy. 

We  have  briefly  investigated  the  solutions  for  the 
sphere  with  k  =1.2  x  10    in  the  neighborhood  of  the  first  few 
eigenvalues.   The  results  are  sketched  in  Figure  11.  The  lowest 
eigenvalue  is  symmetric.   The  solution  branch  which  bifurcates 
from  this  eigenvalue  has  the  same  qualitative  behavior  as  the 
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bifurcation  branch  from  the  lowest  symmetric  eigenvalue  for 
k  =  10   and  k  =  10   .   The  second  eigenvalue  Is  unsymmetrlc 
and  the  solutions  bifurcate  up  and  have  negative  energy  near 
the  eigenvalue.   The  solutions  branch  down  from  the  third 
eigenvalue  P  =  ?-,£-•   They  are  qualitatively  similar  to  the 
bifurcation  branch  from  the  lowest  eigenvalue  for  k  =  10  -^ 
and  k  =  10"^. 


5 . 3   Boundary  Layers . 

We  have  observed  in  Section  5.1,  see  Figures  3  and  ^j 
that  for  the  bifurcation  branch  from  the  lowest  symmetric 
eigenvalue,  a  boundary  layer  forms  at  9  =   Tr/2.   A  second 
boundary  layer  was  observed  to  develop  at  0  =  0  as  k  ->■  0; 
the  solution  is  small  everywhere  except  at  9  =  0 .   In  Figures 
12  and  13  we  show  the  formation  of  this  boundary  layer  by 
presenting  graphs  of  solutions  for  a  sequence  of  decreasing 
values  of  k;  cf .  with  Figures  5  and  4.  These  results  partially 
substantiate  Friedrichs'  conjectures  [4]  about  the  formation 
of  "dimples"  as  k  ->-  0. 

6 .   Remarks  about  Buckling. 

We  shall  discuss  the  relevance  of  the  results  previously 
described  to  the  axisymmetric  buckling  of  perfect  spherical  shells 
under  idealized  experimental  conditions;  that  is  "gedanken 
buckling."   We  discuss  primarily  the  thick  sphere  since  our 
results  are  most  comprehensive  in  this  case.   But  our  remarks  are 
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equally  applicable  to  the  thin  sphere  and  no  doubt  to  a  wide 

range  of  k  values.   Our  basic  premise  is  that  we  have  determined 

all  of  the  relevant  axially  symmetric  equilibrium  states  of  the 

shell  for  all  loads  in  some  large  range :   0  <  P  <  P 

^  ^  —       —     max 

Constant  reference  to  Figures  la  and  2   and  Table  II  is 

required  for  this  discussion. 

We  imagine  the  pressure  increasing  continuously  from 

P  =  0  and  allow  only  axisymmetric  deformations.   The  sphere 

contracts  radially  (corresponding  to  the  trivial  solution  w  =  0) 

for  small  P  since  this  is  the  only  existing  equilibrium  state. 

However,  as  soon  as  P  exceeds  Pf-  -,  ~   '^ f.    n   -    •01045,  two  additional 

equilibrium  states  exist  corresponding  to  the  branches  B^  ,  and 

B^  „.   Suppose  for  the  present  that  the  shell  does  not  jump  or 

is  somehow  constrained  from  jumping  to  another  equilibrium  state. 

Then  as  the  load  continues  to  increase  the  sphere  uniformly 

contracts,  regardless  of  the  numerous  states  which  become  possible, 

until  P  finally  reaches  the  lowest  eigenvalue.   In  the  thick 

sphere  case  this  is  P^-.   At  this  pressure  it  may  be  possible 

for  the  sphere  to  deviate  from  perfect  sphericity  without 

jumping  by  merely  deforming  continuously  into  the  bifurcated 

mode  of  solution.   In  general  if  the  bifurcation  branch  branches 

up  (down)  this  smooth  transition  can  occur  only  for  P  increasing 

(decreasing).   Thus  since  B^-  ,  branches  down  from  P^  such  a 

smooth  transition  cannot  occur  at  P^  as  P  increases. 

5 


IE _c: 

For  other  values  of  k,  e.g.  k  =  1.2x10  ^  the  lowest  eigenvalue 

is  symmetric;  see  Figure  13.   Therefore,  a  smooth  transition 
into  the  bifurcation  branch  as  P  increases  is  possible. 
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We  recall  that  for  the  bifurcation  branches  Bp   ., 
passing  through  each  symmetric  eigenvalue  Pp  ,  the  energy 
Sp   -1  (P)  is  monotone  decreasing.   Thus  as  P  Increases  beyond 
some  Pp   the  unbuckled  state  has  more  energy  than  the  state 
on  Bp   -,  .   If  we  make  the  reasonable  assumption  that  at  a 
bifurcation  point  states  with  least  energy  are  preferred  then 
the  thick  sphere  smoothly  deforms  into  the  state  on  B^  ,  as  P 
goes  from  Pj-  past  P^.   However,  suppose  the  shell  remains 

spherical  as  P  increases  beyond  P^  and  Pk.   Then  at  the 

+ 
unsymmetric  eigenvalue  P„,  whose  branches  B„  .,  bifurcate  up, 

the  smooth  transition  can  occur  onto  one  of  the  unsymmetric 

branches.   Thus  in  summary,  if  jump-buckling  is  somehow  prevented 

or  just  does  not  occur,  the  uniformly  contracted  sphere  smoothly 

deforms  into  a  nonspherical  shape  as  P  increases  past  a  symmetric 

eigenvalue  or  an  unsymmetric  eigenvalue  which  branches  up.   We 

call  smooth  transitions  to  nonspherical  states,  as  P  increases, 

"bifurcation  buckling."   Then  all  the  symmetric  eigenvalues,  Pp  , 

are  bifurcation  buckling  loads  and  so  are  those  unsymmetric 

eigenvalues,  Pp   -,  ,  with  upward  branching  solutions. 

Let  us  now  consider  that  bifurcation  buckling  has 

occurred  at  P  =  P   and  P  continues  to  increase.   The  states  of 

deformation  correspond  to  solutions  on  some  upward  bifurcating 

branch  B   ,.   However,  as  P  increases  beyond  the  upper  critical 

point  P   ^  the  deformation  cannot  vary  continuously  with  P 

since  the  branch  B   -,  ceases  to  exist.   Thus  the  shell  must 

n ,  1 

jump  to  some  other  non-adjacent  equilibrium  state.   We  therefore 
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call  P^  n  a  jump-buckling  load.   Clearly,  V ^^   ^   are  all  Jump- 
buckling  loads  and  so  are  those  P2m+1  1  ^^^^'^  ^^^  endpoints  of 
upward  bifurcating  branches.  For  example,  Py  i  -  0-0858  for 
the  thick  sphere  is  a  jump-buckling  load  while  P   -j^  =  P^ 

=  0.0675  is  not. 

Let  us  return  to  the  unconstrained  thick  sphere  experi- 
ments in  which  bifurcation  buckling  need  not  precede  jump- 
buckling.   As  has  been  observed,  only  uniform  contraction  can 


6,1 


occur  for  P  <  P^  -,  which  we  therefore  term  the  lower  buckling 


load  P^.   Since  P^  is  not  a  bifurcation  buckling  load  we 

5 

conclude  that  as  P  increases  over  the  interval  P   -  0.010^5 
<  P  <  p^  0.07062  the  thick  sphere  can  become  nonspherical 
only  by  jumping  into  a  buckled  state.   The  multiplicity  of 
solutions  varies  over  this  interval.   For  example,   at 
P  =  0.015  there  are  eight  such  solutions  corresponding  to  the 
branches  :  Bg  „;  Bg  ^;  B^  -j_ ;  B~  ^•,    B"  „.   Of  these  solutions 
only  those  on  B^  -,  and  B"  ,  have  more  energy  than  the  unbuckled 
state  at  P  =  O.OI5.   The  states  with  negative  energy  are  in 
some  sense  preferred  to  the  unbuckled  state.   Jumping  to  a  state 
of  negative  energy  cannot  occur  for  V  <  '? r   q.   We  must  invoke 
some  mechanism  for  the  shell  to  jump  into  one  of  these  states 
of  negative  energy.   The  usual  such  mechanism  is  attributed  to 


As  another  example  consider  P  =  0.01996  where  there  are  16 
solutions  on  the  branches:   Bg  „;  Bg  -,  ;  B^  .;  j  =  1,2,..., 7. 

All  but  the  three  solutions  on  B^  ,  and  B±  ,  have  less  energy 

6,1      5,1 

than  the  uniformly  contracted  sphere  with  P  =  0.01996.   In 
Figures  l4  a,b  we  show  ^  of  these  solutions. 
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imperfections  and/or  disturbances  of  the  original  spherical 
shell,  of  the  loading  procedure,  of  the  shell  support,  etc. 
The  transition  is  actually  dynamic  and  since  there  are  many 
buckled  states  present  with  energy  less  than  the  unbuckled 
state  we  do  not  know  which  state  is  the  one  finally  preferred. 
The  number  and  character  of  the  buckled  states  which  coexist 
at  a  given  pressure  may  Influence  the  possibility  of  jump- 
buckling.   Thus  the  triggering  mechanism  which  causes  the  shell 
to  jump,  no  doubt,  depends  on  the  interaction  of  the  imperfec- 
tions, disturbances  and  the  kind  of  buckled  states  that  exist 
at  a  given  pressure. 

The  possible  buckling  mechanism  just  described  is  such 
that  precise  jump-buckling  loads  cannot  be  defined.   However, 
we  note  that  for  P  just  below  a  symmetric  eigenvalue  or  an 
unsymmetric  eigenvalue  whose  bifurcated  solutions  branch  down 
there  must  always  be  buckled  states  very  close  to  the  unbuckled 
state  but  with  slightly  more  energy.   An  example  is  P^  for  the 
thick  sphere.   These  states  could  also  supply  the  triggering 
mechanism  for  jump  buckling  in  experiments  with  very  small 
imperfections  or  disturbances.   In  such  a  case  the  observed 
jump  buckling  load  should  be  less  than  the  lowest  eigenvalue, 
that  is,  P|-  for  the  thick  sphere.   Of  course,  if  the  imperfec- 
tions are  so  small  that  they  do  not  cause  jump  buckling,  then 
as  P  increases  bifurcation  buckling  must  occur  at  P   and  finally 
jump  buckling  must  occur  at  P   ,  . 

We  have  only  considered  axisymmetrlc  deformations  in 
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this  study.   If  we  allow  axlally  unsymmetric  states  there  may 
be  additional  bifurcation  branches  and  other  buckled  states. 
A  numerical  study  of  these  additional  bifurcation  branches  has 
not  yet  been  initiated.   However  an  analysis  of  the  axially 
unsymmetric  bifurcation  branches  near  the  eigenvalues  and  the 
effects  of  imperfections  is  given  in  [8]  using  a  different  shell 
theory . 


-36- 


References 

1.  S.  Timoshenko,  Theory  of  Elastic  Stability.  McGraw-Hill, 
New  York,  1936. 

2.  T.  von  Karman  and  H.  S.  Tsien,  The  Buckling  of  Spherical 
Shells  by  External  Pressure,  J.  Aero.  Sci.  7  (1939), 

pp.  43-50. 

3.  H.  S.  Tsien,  A  Theory  for  the  Buckling  of  Thin  Shells, 
J.  Aero.  Sci.  9  (19^2),  pp.  373-384. 

4.  K.  0.  Friedrichs,  On  the  Minimum  Buckling  Load  for 
Spherical  Shells,  von  Karman  Anniversary  Volume,  Calif. 
Inst.  Tech.,  Pasadena,  1941,  pp.  258-272. 

5.  B.  Altshuler,  Nonlinear  Buckling  of  a  Spherical  Shell, 
Doctoral  thesis.  New  York  Univ.,  1953. 

6.  A.  G.  Gabril'iants  and  V.  I.  Feodos'ev,  Axial ly-Symmetric 
Forms  of  Equilibrium  of  an  Elastic  Spherical  Shell  under 
Uniformly  Distributed  Pressure,  PMM,  25  Cl96l),  pp.  1091- 
1101.  (Engl,  transl.  in  J.  Appl.  Math.  Mech.  25  (1961), 
pp.  1629-1642.) 

7.  T.  Koga  and  N.  J.  Hoff,  The  Axisymmetric  Buckling  of 
Initially  Imperfect  Complete  Spherical  Shells,  Int.  J, 
Solids  and  Structures  5  (1969),  pp.  679-697- 

8.  W.  T.  Koiter,  The  Nonlinear  Buckling  Problem  of  a 
Complete  Spherical  Shell  under  Uniform  External  Pressure, 
Rep.  No.  412,  Lab.  of  Eng.  Mech.,  Technological  Univ., 
Delft,  1968. 


-37- 


9.  J.  M.  T.  Thompson,  The  Rotationally  Symmetric  Branching 
Behaviour  of  a  Complete  Spherical  Shell,  Proc.  Koninkl . 
Nederl.  Akad.  van  Wet.   Amsterdam,  Sec.  B,  6  (1964), 
pp.  295-311. 

10.  E.  L.  Reiss,  H.  J.  Greenberg  and  H.  B.  Keller,  Nonlinear 
Deflections  of  Shallow  Spherical  Shells,  J.  Aero.  Sci. 
24  (1957),  pp.  533-543. 

11.  S.  Timoshenko,  Theory  of  Plates  and  Shells,  McGraw-Hill, 
New  York,  1940. 

12.  E.  Bromberg  and  J.  J.  Stoker,  Nonlinear  Theory  of  Curved 
Elastic  Sheets,  Q.  Appl.  Math.  3  (1945),  pp.  246-265. 

13.  T.  M.  MacRobert,  Spherical  Harmonics,  Dover,  New  York,  1948. 

14.  H.  B.  Keller,  Numerical  Methods  for  Two  Point  Boundary  Value 
Problems ,  Ginn-Blaisdell ,  Waltham,  Mass.,  1968. 

15.  W.  Langford,  Bifurcation  Theory  of  Thin  Spheres, 
in  preparation. 


-38- 


Appendix  I .   Derivation  of  the  Shell  Theory. 

The  two  dimensional  Hooke's  law  for  the  normal  stresses 
Oq   and  a^j,  is 

(A.l)       Oq   =  -^   (Eg  +  VE^)  ,    o^   =   -^  (E^  +  vEg)  . 

1-v  1-v 

Substituting  (2.7)  and  (2.8)  into  (A.l)  yields 

(A. 2a)     ag(r,9)  =  ES(e)  -  z(^)M  ,    ajr,9)  =  ET(e)  -  z(^)N  , 

h  h 


where 


(A. 2b)  , 


and 


/<  S{e)    =    (1-v^)  "^(eg+ve^)  =  (1-V^)  ^[U'-W  +  l  V^+v(U  cot  ©-W)], 
[  T(e)  =  (1-V^)"^(e  +ve0)  =  (l-v^)"^[U  cot  0-W+v{U' -W  +^  V^) ] , 


rn{e)   =  -^  {Kq  +  vfc^)   =  —^   (V  +  vV  cot  0)  , 
(A. 2c)  1-V  1-V 

I  N(9)  =  -^    (k.  +  VK.)  =  -^  (V  cot  0  +  vV)  . 
1-V^  1-V 

The  parameter  k  is  defined  in  (2.9c).   The  dimensionless  quantities 
S  and  T  are  proportional  to  Oq   and  a^  evaluated  on  the  middle 
surface,  z  =  0.   They  are  called  the  membrane  stresses.   The 
dimensionless  quantities  M  and  N,  which  we  call  the  bending  moments, 
are  related  to  the  physical  bending  moments  Mg  and  M.  by 

M  =  -  Mg/(2EhR)  ,     N  =  -  M^/(2EhR)  . 
The  strain  energy  of  the  shell  is  given  approximately  by 
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R+h   2tt 


(A. 3) 


^^1 


(OqEq    +   a.E.)r      sin  6  d6  dcf)  dr  . 


R-h 


0    0 


We  Insert  (2.4),  (2.7),  (2.8),  (A.l)  and  (A. 2)  Into  (A. 3), 
perform  the  r,  (p   integrations  and  use  h/R  <<  1  to  obtain 


TT 


(A. 4)       I    =    2TTR^Eh 


[(Seo+   Te.)    +   R(Mk,+   NK,)]sin   9    dG    . 

o  (p  t)  (p 


(A. 5) 


If  p  is  a  pressure  then  the  work  done,  F,  is  defined  by 


F  =  -pA 


where  A  is  the  change  in  the  volume  enclosed  by  the  mldsurface 
of  the  shell  before  and  after  deformation.   The  pressure  is 
positive  when  it  is  directed  inwards.   To  determine  A  we  first 

denote  the  projections  on  the  equatorial  plane  and  on  the  polar 

* 

axis   of  the  radius  from  the  origin  to  any  point  on  the  shell's 

mid-surface  respectively  by  p  and  6.   Thus  for  the  undeformed 
shell  we  have 

(A. 6)         p  =  p^  =  R  sin  0,     5  =  6^  =  R  cos  9  , 

and  for  the  deformed  shell  we  have 

(A. 7)        p  =  p^  =  R  sin  9  +  H,   6  =  6^  =  R  cos  0  -  J  . 

Here  H  and  J,  the  horizontal  and  vertical  components  of  the 
displacement  of  a  point  on  the  shell's  mid-surface,  are  given  by 


1 

The  polar  axis  is  the  diameter  through  the  north  and  south  poles. 

We  assume  the  shell  is  oriented  so  that  the  polar  axis  is 
vertical.   The  equatorial  plane,  defined  by  9  =  77/2,  0  <_  (j)  <  2tt  , 
is  thus  horizontal. 
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R 


(A. 8)    I  =  CU  cot  e  -W)sln  9  =  e,  sin  9,   i  =  W  cos  9+  U  sin  9. 
J  Is  positive  downwards.   A  is  defined  by 


A     -     T\ 


(pbl    -    Ph^h)  d9  . 


(A. 9) 

We  substitute  (A.6)-(A.8)  into  (A. 9)  and  use  (2.8a)  to  get 


0 


(A. 10)   A  =  ttR- 


[e^  +  3e,+3e^)  +  (U'+W  cot  9 )  (1+e  )  ^]sln^9  df 


0 


Since  |e,|  <<  1,  A  is  given  approximately  by 


■n 


(A. 11)     A  =  ttR- 


[U'  +  W  cot  9  +  3e,]  sin^9  d9  . 


Employing  the  expression  for  e,  given  by  (2.8a)  in  (A. 11),  we 


obtain  after  Integration  by  parts, 

ir 
(A. 12)  A  =  -2uR^  I  W  sin  9  d9   . 

0 

Thus  from  (A.  5)  and  (A. 12)  the  work  done  by  the  normal  pressure 

is 


(A. 13) 


F  =  27tR-^p 


r 


0 


W  sin  e  d9 


We  observe  that  (A. 13)  is  also  the  work  done  by  a  uniformly  . 
distributed  surface  force  p  which  is  directed  toward  r  =  0. 
Hence,  for  this  shell  theory,  normal  pressure  and  a  centrally 
directed  surface  load  have,  within  terms  of  higher  order,  the 
same  work  expression.   Consequently  the  boundary  value  problems 
will  be  the  same  for  both  loads. 
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The  potential  energy  G  is  defined  as 

(A. 14)  G  =  I  -  F  , 

where  I  and  F  are  given  by  (A. 4)  and  (A. 13).   Employing  (2.6), 
(2.8)  and  (A.2b,c)  in  (A.l4),  we  express  G  as  a  functional  of 
u  and  w.   The  Euler  equations  of  this  functional  are 

(A. 15a)      Q  sin  6  +  (S  sin  9)'  -  T  cos  0=0, 

(A. 15b)     (Q  sin  6)'  -  (S+T)sin  9  -  2P  sin  0=0, 

where  the  shear  Q  is  defined  by 

(A. 15c)   Q  =  [(M  sin  0)'  -  N  cos  0  -  SV  sin  0]csc  0, 

and  P  is  defined  in  (2.9c).   Since  the  solutions  are  to  be 
bounded  and  axlsymmetric  with  respect  to  the  polar  axis  the 
boundary  conditions  at  9  =  0,tt  are 

(A. 16)  U  =  W  =  Q  =  0  , 

We   observe   that 

(A. 17)      U=V=M=N=Q=0,      S=T=-P,      W=WqS    (l-v)P    , 

is  a  solution,  for  all  P  and  k,  of  the  system  of  equations 
(A. 2),  (A. 15) J  (2.6)  and  satisfies  the  boundary  conditions 
(A. 16).   It  represents  the  uniformly  contracted  sphere  or 
the  unbuckled  state.   We  seek  other  solutions  in  the  form 
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(A. 18) 


w  =  w  +  wCe) 


s  =  ~p  +  s(e)  ,   T  =  -p  +  t(e)  , 


u  =  u(e)  ,   V  =  v(9)  ,  ...,  Q  =  q(e)  . 


We  substitute  (A.18)  into  (A. 2),  (A. 15),  (2.6)  and  (A.I6)  and 
find  that  the  lower  case  variables  must  satisfy 

(A. 19a)    q  sin  9  +  (s  sin  9)'  -  t  cos  9=0, 

(A. 19b)    (q  sin  9)'  -  (s+t)sin  9=0, 

(A. 19c)    (m  sin  9)'  -  n  cos  9  -  sv  sin  9  -  q  sin  9  +  Pv  sin  9=0, 


(A.19d)    V  =  u  +  w'  , 


(A.19e)   < 


(A.19f)   < 


t  = 


m 


n  = 


( — ^)[u'  -  w  +  ^  v^  +  v(u  cot  9  -  w)]  , 
1-v 

1  12 

( p)[u  cot  9  -  w  +  v(u'  -  w  +  p-  V  )]  , 

1-v 


1-v 

1 


1-v 
(A.19g)    u  =  w'  =  q  =  0  , 


p-  (v'  +  vv  cot  9)  , 
P  (v  cot  9  +  vv' )  , 


=  0,TT   . 


Equations  (A.19a-f)  are  a  system  of  eight  equations  in 
eight  dependent  variables.   We  eliminate  the  four  variables 
s,  n,  u  and  w  from  this  system  and  obtain  the  four  ordinary 
differential  equations  (2.9a,b).   The  boundary  conditions  (2.14) 
follow  from  (A.19d,g).   When  a  solution  of  this  boundary  value 
problem  is  determined  we  evaluate  the  eliminated  variables  from 


(A. 20) 


s  =  q  cot  9  ,    n  =  vm  +  kv  cot 
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The  equations  to  evaluate  u  and  w  are  obtained  by  first 
eliminating  w  from  (A.19e).   This  yields  the  first  order 
ordinary  differential  equation  for  u,  (2.18). 

If  Gp.  Is  the  potential  energy  of  the  unbuckled  state 
(A. 17) J  then  the  relative  potential  energy  e  of  a  buckled 
state  Is 

(A. 21)  e  5  G  -  Gq  . 

We  Insert  (A. 4)  and  (A. 13)  In  (A.l4)  and  express  G  In  terms 
of  m,  q,  t,  V  and  w.   Integration  by  parts  and  CA.19g)  gives 
(2.19). 
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Figure  la 

Graphs  of  P  vs.  the  amplitude  A,  which  Is  defined  In  Equation 
(5.1),  corresponding  to  the  7  lowest  symmetric  eigenvalues  and 
che  4  lowest  unsymmetrlc  eigenvalues  for  the  thick  sphere 
k  =  .001.   The  circles  on  the  P  axis  Indicate  the  symmetric 
eigenvalues  and  the  crosses  the  unsymmetrlc  eigenvalues.   Several 

of  the  lower  symmetric  branches  B   .  and  the  lower  unsymmetrlc 

± 
branches  B   .  are  explicitly  labeled  In  the  graph.   Several  of 

the  upper  P   .  and  lower  P   .  critical  loads  are  also  shown.   The 
n,j  n,j 

branches  B^-  .,  J  =  1,2,..., 7  and  B/-  ,  and  B/-  ^   are  not  completely 
shown  In  the  figure  since  the  corresponding  values  of  A  are  too 
large.   They  are  shown  In  Figure  lb.   The  circles  on  the  graphs 
Indicate  the  pressures  at  which  unsymmetrlc  branches  merge  with 
symmetric  branches. 
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Figure  lb 


Graphs  of  P  vs.  A  for  Bg  ^  with  P  <  Pg,  Bg  ^  and  B   ., 

j  =  1,2,..., 7-   The  scale  of  the  graph  Is  Insufficient  to 

+ 
resolve  the  branches  B^-  .,  j  =  2,  .  .  .  ,7  completely.   The 

+ 
unsymmetric  branches  B^  „  merge  with  the  symmetric  branch 

B/-  ,  at  P  =  PZ"  „  =  0.0124.   The  circle  on  the  graph 
6,1         5,7 

± 
Indicates  the  pressure  at  which  Bj-  ,-,  merge  with  B/-  -,  . 
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Figure  2a 

Graphs  of  (e/2)  ^^  vs.  P  for  the  thick  sphere,  k  =  .001.  The 
P  axis  is  horizontal.  The  circles  on  the  P  axis  Indicate  the 
symmetric  eigenvalues  and  the  crosses  the  unsymmetric  eigen- 
values. Some  of  the  critical  loads  and  branches  are  labeled. 
The  circles  on  the  graphs  indicate  the  pressures  at  which  un- 
symmetric branches  merge  with  symmetric  branches.   The  equal 

■p 
energy  loads  P   .  occur  where  the  curves  cross  the  P  axis  at 
^         n,j 

points  other  than  eigenvalues.   The  P   .  are  too  close  to  P 
^  ^  n,j  n 

to  be  clearly  shown  on  B   .  with  n  >   12.   The  branches  B/-  ^ 

n,j        -  6,0 

and  B(-  .,  J  =  2, 3,..., 7  are  not  completely  shown  in  this 

1/3 
figure  since  the  values  of  -(e/2)  ^^   are  too  large.   They  are 

shown  in  Figure  2b. 
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Figure  2b 

1/3 
Graphs  of  (e/2)  '  -^   vs.  P  for  the  lower  part  of  the  branches 

bifurcating  from  the  lowest  two  eigenvalues  P^-  and  P/-.   The 

eigenvalues  are  not  shown  on  the  graph  since  their  magnitudes 

exceed  the  scale.   The  branch  Ba  „  has  been  continued  as  far 

0,0 

1/^ 
as  P  =  .1,  where  (e/2)  '  ^  ^   -\A.      N(P)  denote  the  number  of 

buckled  solutions  for  P  in  the  interval  P_  <  P  <  P.   The 

following  table  shows  how  N  varies  for  the  range  of  P  given 

in  the  figure. 


i  p 

pL 

ipL 

pL 
5,1 

P^ 
5,3 

'pL 
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pU 

'pU 
/5,2 

'^l6 

p 

0 

p^ 
1  6,1 

i 

pL 
'5,7 

1 

pL 
^5,1 

pL 
5,3 

P^ 
5,5 

■pU 
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0 
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1 
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12 

16 

12 
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Figure  5 

Graphs  of  the  radial  displacement  W(9)  =  w(9)  +  W„,  where  W„  is 
the  uniform  radial  contraction  of  the  unbuckled  state.   The 
solutions  are  for  a  sequence  of  values  of  P  as  indicated  in  the 
following  table.   This  table  also  applies  to  Figure  4. 


Graph  ;  ., 
Number 

1 
2 

3 

K 

5 

J. 

6 

1 

7 

8  i  9 

10  i  11 

12 

13 

P 

^6 

.06 

.04 

.03 

.02 

.015': 

.013 

.011 

.0105 

.0105 

.03 

.075 

.1 

"Rt^  n  nln 

B^ 

, 

=6,1 

■R^ 

J 

> 

^6,0 

t 

I 

^6,0 
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Figure  4 

Graphs  of  the  reduced  circumferential  stress  t(9),  see 
Eq.  (A. 18).   The  solutions  are  for  the  sequence  of  P 
values  given  in  the  table  in  the  caption  for  Figure  3. 
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Figure  5 

Graphs  of  W(9)  for  a  sequence  of  values  of  P  to  show  the  transi- 
tion m  buckled  shapes  from  B^^^  to  Bg^^  ^o  B^  ^.   The  graphs 
numbered  1  and  4  are  for  P  near  the  eigenvalues  Pg  and  P..   The 

graphs  numbered  2  and  5  are  for  P  near  the  critical  loads  p!^ 

T  0,1 

and  P^  ^. 
6,2 


Graph  Number 


Branch 


070   .0712594 


B 


.U 


^  "6,1  =^6,2 


B. 


6,1  I  ^6,1 


or  B 


6,2 


3 


0.678672 


=  P 


6,2 


=  p 


^,1 


.074 


^6,2  °^B^,i  \^^,i 
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Figure  6 

Graphs  of  t(0)  for  the  same  sequence  of  values  of  P  as 
given  in  the  table  of  the  caption  for  Figure  5. 
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Figure  7a 

Graphs  of  W(0)  corresponding  to  points  on  B^  _  for  a 
sequence  of  decreasing  P  values.   Graphs  1,  2  and  3  are 
respectively  for  P  =  .066,  .06l,  .040. 
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Figure  7t) 

Graphs  of  W(0)  for  a  sequence  of  values  of  P  on  B^  .  with  j 

Increasing.   The  purpose  of  the  figure  Is  to  show  how  the  un- 

syiranetric  solutions  gradually  merge  with  the  symmetric  solution 

at  P  =  P^  ^  =  .012397.   Graphs  4,  5,  6  and  7  are  respectively 

for  P  =  .0128125  =  P^  ,  =  P^  ^,  P  =  .0173426  =  P^  ^  =  P^  „, 

bi-J-    D>^  5,3    5,4' 

P  =  .020854  =  P^  ^  =  P^  „,  and  P  =  .013.  The  curve  J   Is  on  b;^  ^ 

Dr'C'  D,  (  ^,  f 

near  P    where  the  unsymmetrlc  and  symmetric  solutions  merge. 
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Flecure  oa 
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Graphs  of  t(9)  corresponding  to  points  on  B,-  -,  for  the 
same  values  of  P  as  In  Figure  Ja.. 
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Figure  8b 

Graphs  of  t(9)  corresponding  to  the  same  solutions  and  pressures 
as  In  Figure  Tb  and  in  addition  Graph  8  for  P  =  .012399'   Since 
P  =  .012599  is  close  to  P^   =  .0123971  this  solution  is  nearly 
symmetric . 
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Figure  9 

Graph  of  P  vs.  A  for  the  solution  branching  from  P-if^-j  the 
lowest  symmetric  eigenvalue  of  the  thin  sphere. 
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Figure  10a 

Graphs  of  (e/2)  '^   vs.  P  for  the  thin  sphere,  k  =  10"^,  for  the 
solutions  branching  from  the  lowest  5  symmetric  eigenvalues. 
As  in  Figure  2a,  the  P  axis  Is  horizontal  and  the  circles  on 
the  P  axis  indicate  the  eigenvalues.   Some  of  the  branches  and 
critical  loads  are  labeled.   The  branch  B, /-  ,  is  not  completely 
shown  in  the  figure.   It  is  shown  in  Figure  10b. 
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FlKure  10b 


Graphs  of  (e/2)  ^^  vs.  P  for  the  lower  part  of  B^g  -,  and 
for  the  branches  B-,/-  „,  B,/-  _^  and  a  part  of  B,/-  _p. 


-78- 


;c/2) 


p 


L 
16, 


\ 


B,a, 


PxlO' 


■J    : 


.6 


-.20 


-.10 


\  .1. 

'   iG,0 


io-^-^ 


-.30.- 


-79- 


Figure  11 

A  sketch  indicating  the  local  behavior  near  the  lowest 
eigenvalues  for  the  sphere  with  k  =  1.2x10  -^ .      The  numeri- 
cal values  for  the  critical  pressures  are  P-^g  ^    .0073217, 
^16,1  =  ^16, 2  =  -00^5567,  P^6^2  ^  -0073552,  P^^  =  .007^640, 
P^^  =  .007392. 
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Figure  12 

Graphs  of  W(0)  for  P  =  6  x 10   and  a  sequence  of  values  of 
k  decreasing  from  k  =  10"-^  (the  thin  sphere)  to  k  =  2  xlO~  . 
For  k  =  10   the  solution  corresponds  to  a  point  on  B^ /-  -,. 
The  values  of  kxlO-^  are  1,  .8,  .6,    A,    and  .2  for  the  graphs 
numbered  1  through  5  in  that  order. 
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Figure  1^ 

Graphs  of  t(9)  for  P  =  6x10"   and  the  same  sequence  of 
values  of  k  as  in  Figure  12. 
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Figure  l4a 

Graphs  of  W(e)  for  two  of  the  l6  solutions  at  P  =  .OI996 
for  the  thick  sphere.   Two  other  solutions  are  shown  in 
Figure  l4b.   The  remaining  solutions  are  not  shown  since 
they  differ  little  from  the  first  solution  shown  in 
Figure  l4b. 


-86- 


-87- 


Figure  l4b 
Graphs  of  W(9)  for  two  of  the  l6  solutions  at  P  =  .OI996. 
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^A  nonlinear  thin  shell  theory  is  derived  for  the  axisymmetric 
buckling  of  spherical  shells  subjected  to  either  a  pressure  or  a 
centrally  directed  surface  load.   The  theory  is  reduced  to  a  boundary 
value  problem  for  a  system  of  four  first  order  ordinary  differential 
equations.   Numerical  solutions  of  this  boundary  value  problem  are 
obtained  by  the  shooting  and  parallel  shooting  methods.   An  extensive 
numerical  study  is  made  of  the  nonlinear  deformations  of  the  shells. 
We  find  for  example,  that  all  solution  branches  that  bifurcate  from 
the  eigenvalues  of  the  linearized  buckling  theory  are  connected  to 
each  other  by  means  of  intermediate  branches.   Some  implications  of 
the  numerical  results  concerning  the  buckling  of  spherical  shells  are 
discussed. 
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